library(BMA)
library(foreign)
library(BAS)

## Run state do file first
abddata<-read.dta("immigration_data.dta")

## G Prior
abd.gprior = bas.lm(h1bvisas_scale~groupcat2+groupcat1+dscore2+female+age+married+educ+racewhite+income2+pid+techzip2,data=abddata, n.models = 2^11, alpha=4, prior="g-prior", initprobs="Uniform", update=NULL, bestmodel = NULL, bestmarg = NULL)

abd.gprior
plot(abd.gprior)
summary(abd.gprior)
image(abd.gprior, subset=-1, vlas=0)
abdgprior.coef = coefficients(abd.gprior)
abdgprior.coef
plot(abdgprior.coef)

## Hyper-g
abd.hypergprior = bas.lm(h1bvisas_scale~groupcat2+groupcat1+dscore2+female+age+married+educ+racewhite+income2+pid+techzip2,data=abddata, n.models = 2^11, alpha=4, prior="hyper-g", initprobs="Uniform", update=NULL, bestmodel = NULL, bestmarg = NULL)

abd.hypergprior
plot(abd.hypergprior)
summary(abd.hypergprior)
image(abd.hypergprior, subset=-1, vlas=0)
abdhyperg.coef = coefficients(abd.hypergprior)
abdhyperg.coef
plot(abdhyperg.coef)

## Hyper-g Laplace
abd.hyperglprior = bas.lm(h1bvisas_scale~groupcat2+groupcat1+dscore2+female+age+married+educ+racewhite+income2+pid+techzip2,data=abddata, n.models = 2^11, alpha=4, prior="hyper-g-laplace", initprobs="Uniform", update=NULL, bestmodel = NULL, bestmarg = NULL)

abd.hyperglprior
plot(abd.hyperglprior)
summary(abd.hyperglprior)
image(abd.hyperglprior, subset=-1, vlas=0)
abdhypergl.coef = coefficients(abd.hyperglprior)
abdhypergl.coef
plot(abdhypergl.coef)

## ZS
abd.zsprior = bas.lm(h1bvisas_scale~groupcat2+groupcat1+dscore2+female+age+married+educ+racewhite+income2+pid+techzip2,data=abddata, n.models = 2^11, alpha=NULL, prior="ZS-null", initprobs="Uniform", update=NULL, bestmodel = NULL, bestmarg = NULL)

abd.zsprior
plot(abd.zsprior)
summary(abd.zsprior)
image(abd.zsprior, subset=-1, vlas=0)
abdzs.coef = coefficients(abd.zsprior)
abdzs.coef
plot(abdzs.coef)

## AIC
abd.aicprior = bas.lm(h1bvisas_scale~groupcat2+groupcat1+dscore2+female+age+married+educ+racewhite+income2+pid+techzip2,data=abddata, n.models = 2^11, alpha=NULL, prior="AIC", initprobs="Uniform", update=NULL, bestmodel = NULL, bestmarg = NULL)

abd.aicprior
plot(abd.aicprior)
summary(abd.aicprior)
image(abd.aicprior, subset=-1, vlas=0)
abdaic.coef = coefficients(abd.aicprior)
abdaic.coef
plot(abdaic.coef)

## BIC
abd.bicprior = bas.lm(h1bvisas_scale~groupcat2+groupcat1+dscore2+female+age+married+educ+racewhite+income2+pid+techzip2,data=abddata, n.models = 2^11, alpha=NULL, prior="BIC", initprobs="Uniform", update=NULL, bestmodel = NULL, bestmarg = NULL)

abd.bicprior
plot(abd.bicprior)
summary(abd.bicprior)
image(abd.bicprior, subset=-1, vlas=0)
abdbic.coef = coefficients(abd.bicprior)
abdbic.coef
plot(abdbic.coef)

## EB
abd.ebprior = bas.lm(h1bvisas_scale~groupcat2+groupcat1+dscore2+female+age+married+educ+racewhite+income2+pid+techzip2,data=abddata, n.models = 2^11, alpha=NULL, prior="EB-local", initprobs="Uniform", update=NULL, bestmodel = NULL, bestmarg = NULL)

abd.ebprior
plot(abd.ebprior)
summary(abd.ebprior)
image(abd.ebprior, subset=-1, vlas=0)
abdeb.coef = coefficients(abd.ebprior)
abdeb.coef
plot(abdeb.coef)

